function f= twitch(peak_force, rise_time_ms, start_time_ms,end_ms)


t=0:1:999;
f=zeros(1,end_ms+1);
f( 1+start_time_ms:(start_time_ms+1000) )= peak_force*t./rise_time_ms.*exp(1-(t./rise_time_ms));

